Introduction¶
In this project we are going to analyze 'Brain tumor dataset', a T1 weighted MRI dataset from figshare,
extract features and create image classification model based on two Convolutional neural network models that predicts both brain tumor type and it's location.
Contents:
- Data exploration
1.1 - Image extraction
1.2 - Splitting Data
1.3 - Preview images
1.4 - Black/White balance
1.5 - Skull size
1.6 - Class weights
1.7 - Feature scaling - 1st Model (tumor type)
2.1 - Building a model
2.2 - Visualising results
2.3 - Testing the model - Predicting tumor location
3.1 - Mask extraction
3.2 - Splitting Data
3.3 - Preview the masks
3.4 - Tumor size
3.5 - Feature scaling - 2nd Model (tumor location)
- Conclusion
Highly reccomended to read the data-description.txt file to understand data organization.
*The code was written in Python 3.12 later changed to 3.10 due to usage of tensorflow-gpu.
I Data exploration¶
1.1 - Image extraction¶
The Dataset is in .mat format, so we extract the image values to regular image files.
For this task we're going to use h2py and cv2 libraries and iterate through mat dataset.
*We don't really need to do that, but I want to show the process as working with regular images.
Technically we supposed to to assign labels to integer values because later on we'll use loss='sparse_categorical_crossentropy' which is prepared for said value type.
Since we're going to create two different dataset splits, keep in mind that 'train', 'val' and 'test' datasets are for 1st model
while 'train_2', 'val_2' and 'test_2' are prepared for 2nd model.
1st split is by class (tumor type) and second one is by image and mask (tumor location).
import tensorflow as tf
import numpy as np
import os
# set up gpus
gpus = tf.config.experimental.list_physical_devices('GPU')
for gpu in gpus:
tf.config.experimental.set_memory_growth(gpu, True)
# print(gpus)
# tf.device('/GPU:0')
# apply stable seed to avoid inconsistent results every time we test new setup
tf.config.experimental.enable_op_determinism()
tf.keras.utils.set_random_seed(0)
tf.random.set_seed(0)
np.random.seed(0)
# define paths for datasets
mat_path = 'D:\\data-science\\brain-tumor-data\\mat'
data_path = 'D:\\data-science\\brain-tumor-data\\images'
data_path_split = data_path + '_split'
mask_path = 'D:\\data-science\\brain-tumor-data\\mask'
path_to_images = os.path.join(mask_path, 'x-image')
path_to_masks = os.path.join(mask_path, 'y-mask')
matrices = os.listdir(mat_path)
import matplotlib.pyplot as plt
# DON'T MIND THIS, it's dark theme for matplotlib
WHITE_MID = '#b5b5b5'
GREY_DARK = '#141414'
plt.rcParams['figure.facecolor'] = GREY_DARK
plt.rcParams['text.color'] = WHITE_MID
plt.rcParams['axes.facecolor'] = GREY_DARK
plt.rcParams['axes.edgecolor'] = WHITE_MID
plt.rcParams['axes.labelcolor'] = WHITE_MID
plt.rcParams['axes.titlecolor'] = WHITE_MID
plt.rcParams['grid.color'] = WHITE_MID
plt.rcParams['grid.linestyle'] = '--'
plt.rcParams['grid.linewidth'] = 0.5
plt.rcParams['axes.grid'] = True
plt.rcParams['axes.linewidth'] = 1
plt.rcParams['xtick.color'] = WHITE_MID
plt.rcParams['ytick.color'] = WHITE_MID
plt.rcParams['legend.edgecolor'] = WHITE_MID
plt.rcParams['legend.labelcolor'] = WHITE_MID
1.2 - Train, val and test split¶
For this task we use splitfolders library. This makes us sure that data is correctly split without any leaks.
Ratio:
train: 70%
val: 20%
test: 10%
train = tf.keras.utils.image_dataset_from_directory(os.path.join(data_path_split, 'train'), batch_size=8, color_mode='grayscale', image_size=(256, 256), label_mode='categorical', shuffle=True, seed=42)
test = tf.keras.utils.image_dataset_from_directory(os.path.join(data_path_split, 'test'), batch_size=8, color_mode='grayscale', image_size=(256, 256), label_mode='categorical', shuffle=True, seed=42)
val = tf.keras.utils.image_dataset_from_directory(os.path.join(data_path_split, 'val'), batch_size=8, color_mode='grayscale', image_size=(256, 256), label_mode='categorical', shuffle=True, seed=42)
data_iterator = test.as_numpy_iterator()
batch = data_iterator.next()
labels = test.class_names
labels
Found 2144 files belonging to 3 classes. Found 308 files belonging to 3 classes. Found 612 files belonging to 3 classes.
['glioma', 'meningioma', 'pituitary']
To ensure that there are no leaks that might overfit our model learning process or lead to false accuracy results during testing
I applied custom function that takes every image from one split and compares it with each image from another split.
The function applies broadcasting. I tested this by copying images between folders by hand. The resulted as predicted (number of duplicates increased for tested directories).
train_images = tf.concat([img for img, label in train], axis=0)
val_images = tf.concat([img for img, label in val], axis=0)
test_images = tf.concat([img for img, label in test], axis=0)
def compare_imgs(data_1, data_2, chunk_size=10):
nduplicates = 0
for i in range(0, data_1.shape[0], chunk_size):
data_1_chunk = data_1[i : i + chunk_size]
# broadcasting
data_1_chunk_exp = tf.expand_dims(data_1_chunk, axis=1)
data_2_exp = tf.expand_dims(data_2, axis=0)
# Compare
equality_matrix = tf.reduce_all(data_1_chunk_exp == data_2_exp, axis=[2, 3])
# Count duplicates in the chunk
nduplicates += tf.reduce_sum(tf.cast(tf.reduce_any(equality_matrix, axis=-1), tf.int32))
return nduplicates
count = compare_imgs(train_images, val_images)
print(f"There are {count} data leaks between train and val splits.")
count = compare_imgs(train_images, test_images)
print(f"There are {count} data leaks between train and test splits.")
count = compare_imgs(val_images, test_images)
print(f"There are {count} data leaks between val and test splits.")
count = compare_imgs(test_images, test_images)
print(f"There are {count} data leaks between test and test splits.")
There are 0 data leaks between train and val splits. There are 0 data leaks between train and test splits. There are 0 data leaks between val and test splits. There are 308 data leaks between test and test splits.
1.3 - Preview the images¶
Data iterator allows us to move through data and preview it more swiftly.
import seaborn as sb
fig, axes = plt.subplots(nrows=1, ncols=8, figsize=(24, 6))
for i, img in enumerate(batch[0][:8]):
class_name = labels[np.argmax(batch[-1][i])]
axes[i].imshow(img, cmap='gray')
axes[i].set_title(class_name)
axes[i].set_axis_off()
for key, value in zip(['batches', 'width', 'height', 'channels'], batch[0].shape):
print(f"{key:>10}: {value}")
plt.tight_layout()
plt.show()
batches: 8
width: 256
height: 256
channels: 1
Our dataset contains 512px x 512px images and each pixel has 3 channels so basically each image contains 786 432 values,
with color_mode = 'grayscale' we reduce it to 262 144. Also keras resized our files to size of 256px x 256px so we end up with 65 536 values per image.
1.4 - Black/white balance¶
Every point on the plot is a single image, we can note that most of images hang around 70.
Note that we skip black background. We just want to check mean color of head area.
from scipy.stats.mstats import hmean
# skipping values to 10 because it's background and it will lower our mean
all_images = np.concatenate([train_images, val_images, test_images], axis=0)
all_means = [np.mean(x[x > 10]) for x in all_images]
print(f'{"Mean":>14}: {np.mean(all_means)}\n'
f'{"Harmonic mean":>14}: {hmean(all_means)}')
fig, axes = plt.subplots(nrows=2, ncols=1, figsize=(20, 10))
sb.scatterplot(y=[x for x in range(len(all_means))], x=all_means, color=WHITE_MID, s=20, ax=axes[0], marker='.', linewidth=0)
sb.kdeplot(x=all_means, ax=axes[1], color=WHITE_MID)
padding = (min(all_means) - 5, max(all_means) + 5)
axes[0].set_xlim(padding)
axes[0].set_title('All images')
axes[1].set_xlim(padding)
axes[1].set_title('Density')
plt.tight_layout()
Mean: 70.33614349365234 Harmonic mean: 68.72393664459248
1.5 Head area¶
By calculating amount of white pixels inside image we can extract images with specific skull/head area.
threshold is max value that we consider/define as black background. The more the noisy our background is the higher the threshold we apply.
In summary this plot returns images where head fills specified % (see percentage variable) of area in image.
threshold = 10
percentage = 25
head_pixels = []
# sub plot created for 18 images, I checked the amount before
fig, axes = plt.subplots(nrows=3, ncols=6, figsize=(24, 12))
axes = axes.flatten()
small_head = 0
i = 0
for img in all_images:
n_white_pixels = np.sum(img > threshold)
percentage_white_pixels = (n_white_pixels / (256*256)) * 100
head_pixels.append(percentage_white_pixels)
if percentage_white_pixels < percentage:
if i < 18:
small_head += 1
axes[i].set_title(f"white pixels: {round(percentage_white_pixels, 2)} %")
axes[i].set_axis_off()
axes[i].imshow(img, cmap='gray')
i += 1
print(f"There are {small_head} image/s with head pixels filling less than {percentage}%.")
plt.tight_layout()
There are 18 image/s with head pixels filling less than 25%.
This allow us to categorize images based on skull area.
plt.figure(figsize=(16, 6))
sb.histplot(data=head_pixels, color=WHITE_MID, linewidth=1, edgecolor=GREY_DARK)
plt.title('Head area distribution')
plt.xlabel('% of white (head) pixels in image')
plt.show()
1.6 - Class weights¶
files_per_label = lambda folder: os.listdir(os.path.join(data_path, folder))
# generate dictionary with lengths for each label
labels_len = {label: len(files_per_label(label)) for i, label in enumerate(labels)}
sb.barplot(data=labels_len, x=labels_len.keys(), y=labels_len.values(), color=WHITE_MID, linewidth=0)
print(f"{'Class':<10} Count")
for label, count in labels_len.items():
print(f"{label:<10} : {count}")
plt.title('Number of files per class')
plt.show()
Class Count glioma : 1426 meningioma : 708 pituitary : 930
To handle class imbalance, we can calculate weights with given formula:
$w_c = \frac{N}{N_c * C}$
where:
$ w_c $ = class weight
$ N $ = number of samples in dataset
$ N_c $ = number of samples in class
$ C $ = number of unique classes in dataset
class_weights = dict()
for i, count in enumerate(labels_len.values()):
class_weights[i] = round(sum(labels_len.values()) / (len(labels) * count), 1)
print(f"{'Class':<10} Weight")
for i, weight in class_weights.items():
print(f"{labels[i]:<10} : {weight}")
Class Weight glioma : 0.7 meningioma : 1.4 pituitary : 1.1
1.7 - Feature scaling¶
Usually when we work with Images, the max value is 255 so we're going to transform our values from [min, max] = [0, 255] to [0, 1].
print(batch[0].min(), batch[0].max(), sep='\n')
0.0 249.75
train = train.map(lambda x, y: (x/255, y))
test = test.map(lambda x, y: (x/255, y))
val = val.map(lambda x, y: (x/255, y))
iterator = test.as_numpy_iterator()
batch = iterator.next()
min_value = iterator.next()[0].min()
max_value = iterator.next()[0].max()
print(min_value, max_value, sep='\n')
0.0 0.99313724
II 1st Model¶
2.1 - Building a neural network¶
In short our model base on convolutional layers, where each layer contain filter that convolve with input matrix (image) resulting filtered output.
As input shape we apply values that corresponds to:
- width = 256
- height = 256
- channel = 1 (grayscale).
ReLU activation function cleans our output from all unnecessary negative values, replacing them with zeros.
Dropout layers drops random neurons each time convolution layer is being applied, reducing chance to overfit.
from tensorflow.keras.models import Sequential
from tensorflow.keras.layers import Conv2D, MaxPooling2D, Dense, Flatten, Dropout, Input, ReLU
from tensorflow.keras.regularizers import l2
from tensorflow.keras.optimizers import Adam
from tensorflow.keras.callbacks import ModelCheckpoint, EarlyStopping
def create_model(kernel_size=3, learning_rate=0, regulate=0, slope=0):
model = Sequential()
model.add(Input(shape=(256, 256, 1)))
for num_filters in [16, 32, 16]:
model.add(Conv2D(num_filters, (kernel_size, kernel_size), stride=1, activation='relu'))
model.add(Dropout(0.3))
model.add(MaxPooling2D())
model.add(Flatten())
model.add(Dense(256, activation='relu'))
model.add(Dropout(0.5))
model.add(Dense(3, activation='softmax'))
model.compile(optimizer=Adam(learning_rate=learning_rate), loss='categorical_crossentropy', metrics=['categorical_accuracy'])
return model
model = create_model(learning_rate=0.0005)
model.summary()
Model: "sequential_15"
_________________________________________________________________
Layer (type) Output Shape Param #
=================================================================
conv2d_43 (Conv2D) (None, 254, 254, 16) 160
dropout_51 (Dropout) (None, 254, 254, 16) 0
max_pooling2d_42 (MaxPoolin (None, 127, 127, 16) 0
g2D)
conv2d_44 (Conv2D) (None, 125, 125, 32) 4640
dropout_52 (Dropout) (None, 125, 125, 32) 0
max_pooling2d_43 (MaxPoolin (None, 62, 62, 32) 0
g2D)
conv2d_45 (Conv2D) (None, 60, 60, 16) 4624
dropout_53 (Dropout) (None, 60, 60, 16) 0
max_pooling2d_44 (MaxPoolin (None, 30, 30, 16) 0
g2D)
flatten_14 (Flatten) (None, 14400) 0
dense_28 (Dense) (None, 256) 3686656
dropout_54 (Dropout) (None, 256) 0
dense_29 (Dense) (None, 3) 771
=================================================================
Total params: 3,696,851
Trainable params: 3,696,851
Non-trainable params: 0
_________________________________________________________________
n_epochs = 20
model_checkpoint = ModelCheckpoint(
'best_model_1.h5',
monitor='val_categorical_accuracy',
save_best_only=True,
restore_best_weights=True,
verbose=0
)
hist = model.fit(train, validation_data=val, epochs=n_epochs, class_weight=class_weights, callbacks=[model_checkpoint])
Epoch 1/20 268/268 [==============================] - 7s 22ms/step - loss: 0.6757 - categorical_accuracy: 0.6604 - val_loss: 0.5687 - val_categorical_accuracy: 0.7337 Epoch 2/20 268/268 [==============================] - 6s 21ms/step - loss: 0.4001 - categorical_accuracy: 0.8223 - val_loss: 0.4343 - val_categorical_accuracy: 0.8824 Epoch 3/20 268/268 [==============================] - 6s 21ms/step - loss: 0.2592 - categorical_accuracy: 0.8909 - val_loss: 0.3510 - val_categorical_accuracy: 0.8954 Epoch 4/20 268/268 [==============================] - 6s 21ms/step - loss: 0.1895 - categorical_accuracy: 0.9123 - val_loss: 0.2822 - val_categorical_accuracy: 0.9085 Epoch 5/20 268/268 [==============================] - 6s 21ms/step - loss: 0.1508 - categorical_accuracy: 0.9394 - val_loss: 0.2891 - val_categorical_accuracy: 0.9216 Epoch 6/20 268/268 [==============================] - 6s 21ms/step - loss: 0.1190 - categorical_accuracy: 0.9534 - val_loss: 0.2571 - val_categorical_accuracy: 0.9101 Epoch 7/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0900 - categorical_accuracy: 0.9660 - val_loss: 0.2471 - val_categorical_accuracy: 0.9167 Epoch 8/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0832 - categorical_accuracy: 0.9683 - val_loss: 0.2223 - val_categorical_accuracy: 0.9265 Epoch 9/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0709 - categorical_accuracy: 0.9711 - val_loss: 0.1985 - val_categorical_accuracy: 0.9297 Epoch 10/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0537 - categorical_accuracy: 0.9776 - val_loss: 0.1863 - val_categorical_accuracy: 0.9379 Epoch 11/20 268/268 [==============================] - 6s 22ms/step - loss: 0.0430 - categorical_accuracy: 0.9841 - val_loss: 0.1814 - val_categorical_accuracy: 0.9395 Epoch 12/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0306 - categorical_accuracy: 0.9902 - val_loss: 0.2058 - val_categorical_accuracy: 0.9265 Epoch 13/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0226 - categorical_accuracy: 0.9930 - val_loss: 0.1961 - val_categorical_accuracy: 0.9379 Epoch 14/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0575 - categorical_accuracy: 0.9813 - val_loss: 0.2218 - val_categorical_accuracy: 0.9297 Epoch 15/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0400 - categorical_accuracy: 0.9865 - val_loss: 0.1927 - val_categorical_accuracy: 0.9281 Epoch 16/20 268/268 [==============================] - 6s 22ms/step - loss: 0.0321 - categorical_accuracy: 0.9893 - val_loss: 0.1844 - val_categorical_accuracy: 0.9330 Epoch 17/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0175 - categorical_accuracy: 0.9944 - val_loss: 0.1647 - val_categorical_accuracy: 0.9461 Epoch 18/20 268/268 [==============================] - 6s 22ms/step - loss: 0.0200 - categorical_accuracy: 0.9921 - val_loss: 0.2293 - val_categorical_accuracy: 0.9134 Epoch 19/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0322 - categorical_accuracy: 0.9888 - val_loss: 0.2139 - val_categorical_accuracy: 0.9428 Epoch 20/20 268/268 [==============================] - 6s 21ms/step - loss: 0.0161 - categorical_accuracy: 0.9935 - val_loss: 0.1869 - val_categorical_accuracy: 0.9330
2.2 - Visualising results¶
I split it to 4 different plots where green line is the best epoch.
acc_labels = list(hist.history.keys())
best_epoch = hist.history['val_loss'].index(min(hist.history['val_loss']))
print(hist.history['val_loss'][best_epoch],
hist.history['val_categorical_accuracy'][best_epoch],
sep='\n'
)
fig, axes = plt.subplots(nrows=2, ncols=2, figsize=(20, 10))
axes = axes.flatten()
sb.lineplot(data=hist.history, y='categorical_accuracy', x=[x for x in range(1, n_epochs+1, 1)], label='train', ax=axes[0], color=WHITE_MID)
sb.lineplot(data=hist.history, y='val_categorical_accuracy', x=[x for x in range(1, n_epochs+1, 1)], label='validation', ax=axes[0], linestyle='--', color=WHITE_MID)
sb.lineplot(data=hist.history, y='loss', x=[x for x in range(1, n_epochs+1, 1)], label='train', ax=axes[1], color=WHITE_MID)
sb.lineplot(data=hist.history, y='val_loss', x=[x for x in range(1, n_epochs+1, 1)], label='validation', ax=axes[1], linestyle='--', color=WHITE_MID)
sb.lineplot(data=hist.history, y='categorical_accuracy', x=[x for x in range(1, n_epochs+1, 1)], label='accuracy', ax=axes[2], color=WHITE_MID)
sb.lineplot(data=hist.history, y='loss', x=[x for x in range(1, n_epochs+1, 1)], label='loss', ax=axes[2], linestyle='--', color=WHITE_MID)
sb.lineplot(data=hist.history, y='val_categorical_accuracy', x=[x for x in range(1, n_epochs+1, 1)], label='accuracy', ax=axes[3], color=WHITE_MID)
sb.lineplot(data=hist.history, y='val_loss', x=[x for x in range(1, n_epochs+1, 1)], label='loss', ax=axes[3], linestyle='--', color=WHITE_MID)
for i in range(4):
axes[i].legend(loc='lower right', bbox_to_anchor=(1.005, 1), ncol=2)
axes[i].set_xlabel('Epoch')
axes[i].set_ylim([0, 1])
axes[i].set_xlim([1, n_epochs])
axes[i].set_ylabel('')
axes[i].axvline(x=best_epoch+1, color='darkgreen')
axes[0].set_title('Accuracy (train / val)', weight='bold')
axes[1].set_title('Loss (train / val)', weight='bold')
axes[2].set_title('Train (accuracy / loss)', weight='bold')
axes[3].set_title('Validation (accuracy / loss)', weight='bold')
plt.tight_layout()
0.16467131674289703 0.9460784196853638
2.3 - Testing the model¶
Both Training and Validation curves looks very promising, however we have to test our model for final results.
y_true = []
y_pred = []
for batch_images, batch_labels in test:
y_true_batch = np.argmax(batch_labels, axis=1)
y_pred_batch = np.argmax(model.predict(batch_images, verbose=0), axis=1)
y_true.extend(y_true_batch)
y_pred.extend(y_pred_batch)
y_true = np.array(y_true)
y_pred = np.array(y_pred)
from sklearn.metrics import classification_report, confusion_matrix
report = classification_report(y_true, y_pred, target_names=labels)
print(report)
cm = confusion_matrix(y_true=y_true, y_pred=y_pred)
fig = plt.figure(figsize=(4, 4))
fig = sb.heatmap(cm, annot=True, fmt='d', cbar=False, cmap='Greys_r', linewidths=0)
plt.xlabel('Predicted')
plt.ylabel('True')
plt.axis('off')
plt.tight_layout()
precision recall f1-score support
glioma 0.93 0.97 0.95 143
meningioma 0.92 0.83 0.88 72
pituitary 0.97 0.97 0.97 93
accuracy 0.94 308
macro avg 0.94 0.92 0.93 308
weighted avg 0.94 0.94 0.94 308
III Finding tumor location¶
First we prepare the file paths to our images and masks that we extracted before.
print(path_to_images, path_to_masks, sep='\n')
D:\data-science\brain-tumor-data\mask\x-image D:\data-science\brain-tumor-data\mask\y-mask
input_img_paths = sorted(
[
os.path.join(path_to_images, fname)
for fname in os.listdir(path_to_images)
if fname.endswith(".png")
]
)
mask_img_paths = sorted(
[
os.path.join(path_to_masks, fname)
for fname in os.listdir(path_to_masks)
if fname.endswith(".png") and not fname.startswith(".")
]
)
Checking name pairing with shuffled paths for extra accuracy.
import random
random.Random(0).shuffle(input_img_paths)
random.Random(0).shuffle(mask_img_paths)
count = 0
for nth_img_path, nth_mask_path in zip(input_img_paths, mask_img_paths):
img_id = os.path.basename(nth_img_path)[:4]
mask_id = os.path.basename(nth_mask_path)[:4]
if img_id != mask_id:
count += 1
print(f"Misaligned file names: {count}")
len(input_img_paths) == len(mask_img_paths)
Misaligned file names: 0
True
Function for generating dataset.
def get_dataset(
batch_size,
img_size,
X_files,
y_files,
):
"""Returns a TF Dataset."""
def load_img_masks(nth_img_path, nth_mask_path):
input_img = tf.io.read_file(nth_img_path)
input_img = tf.io.decode_png(input_img, channels=1)
input_img = tf.image.resize(input_img, img_size)
input_img = tf.image.convert_image_dtype(input_img, "float32")
mask_img = tf.io.read_file(nth_mask_path)
mask_img = tf.io.decode_png(mask_img, channels=1)
mask_img = tf.image.resize(mask_img, img_size, method="nearest")
mask_img = tf.image.convert_image_dtype(mask_img, "uint8")
return input_img, mask_img
dataset = tf.data.Dataset.from_tensor_slices((X_files, y_files))
dataset = dataset.map(load_img_masks, num_parallel_calls=tf.data.AUTOTUNE)
return dataset.batch(batch_size)
train_nsamples = round((3064 * 0.7))
val_nsamples = round((3064 * 0.2))
test_nsamples = round((3064 * 0.1))
print(train_nsamples, val_nsamples, test_nsamples)
sum([train_nsamples, val_nsamples, test_nsamples]) == len(input_img_paths)
2145 613 306
True
3.2 Splitting Data¶
For X we apply original images and for y, the masks.
X_train_files = input_img_paths[:train_nsamples]
y_train_files = mask_img_paths[:train_nsamples]
X_val_files = input_img_paths[train_nsamples:(train_nsamples + val_nsamples)]
y_val_files = mask_img_paths[train_nsamples:(train_nsamples + val_nsamples)]
X_test_files = input_img_paths[-test_nsamples:]
y_test_files = mask_img_paths[-test_nsamples:]
train_2 = get_dataset(
8,
(256, 256),
X_train_files,
y_train_files
)
val_2 = get_dataset(
8,
(256, 256),
X_val_files,
y_val_files
)
test_2 = get_dataset(
8,
(256, 256),
X_test_files,
y_test_files
)
# checks if the files inside split are unique (no duplicates)
train_files_ids = set(X_train_files)
print(len(train_files_ids) == train_nsamples)
val_files_ids = set(X_val_files)
print(len(val_files_ids) == val_nsamples)
test_files_ids = set(X_test_files)
print(len(test_files_ids) == test_nsamples)
# this should return 3x False for proving there's no leaks between splits
print(bool(train_files_ids.intersection(val_files_ids)))
print(bool(train_files_ids.intersection(test_files_ids)))
print(bool(test_files_ids.intersection(val_files_ids)))
True True True False False False
3.3 - Preview the masks¶
nth_batch = random.randint(0, 10)
fig, axes = plt.subplots(2, 8, figsize=(32, 8))
axes = axes.flatten()
axes[0].set_ylabel('image', size=20)
axes[8].set_ylabel('mask', size=20)
for img, mask in test_2.take(nth_batch):
for i in range(8):
axes[i].imshow(img[i], cmap='grey')
axes[i].set_xticklabels([])
axes[i].set_yticklabels([])
axes[i+8].imshow(mask[i], cmap='grey')
axes[i+8].set_xticklabels([])
axes[i+8].set_yticklabels([])
plt.tight_layout()
3.4 Tumor area¶
This step is similar to head area step. But this time we calculate the amount of white pixels inside mask images.
import cv2 as cv
y_all = np.concatenate((y_train_files, y_val_files, y_test_files))
tumor_pixels = []
for nth_mask_path in y_all:
mask = cv.imread(nth_mask_path)
n_white_pixels = np.sum(mask == 255)
percentage_white_pixels = (n_white_pixels / (512*512)) * 100
tumor_pixels.append(percentage_white_pixels)
plt.figure(figsize=(16, 6))
sb.histplot(data=tumor_pixels, color=WHITE_MID, linewidth=1, edgecolor=GREY_DARK)
plt.title('Tumor area')
plt.xlabel('% of white pixels (tumor) filling mask image')
plt.tight_layout()
3.5 - Feature scaling¶
batch_sample = train_2.take(1)
for img, mask in batch_sample:
print(
f"min: {img[0].numpy().min()},",
f"max: {img[0].numpy().max()}",
)
min: 0.0, max: 215.75
train_2 = train_2.map(lambda x, y: (x/255, y/255))
val_2 = val_2.map(lambda x, y: (x/255, y/255))
test_2 = test_2.map(lambda x, y: (x/255, y/255))
batch_sample = train_2.take(1)
for img, mask in batch_sample:
print(
f"min: {img[0].numpy().min()},",
f"max: {img[0].numpy().max()}",
)
min: 0.0, max: 0.8460784554481506
IV 2nd Model¶
4.1 - Building neural network¶
In first model we used standard convolutional neural network. This time we're also going to use CNN but with U-Net arichtecture.
Also, instead of using Sequential API, we use Functional API, which allows the model to share its layers between the encoder and decoder.
from keras import layers, Model
from keras.metrics import BinaryAccuracy
def create_unet_model(learning_rate=0.0001):
# INPUT
inputs = Input(shape=(256, 256, 1))
x = layers.Conv2D(16, 3, strides=2, padding="same")(inputs)
x = layers.BatchNormalization()(x)
x = layers.Activation("relu")(x)
# ENCODER
for filters in [32, 64, 128]:
x = layers.SeparableConv2D(filters, 3, padding="same", activation='relu')(x)
x = layers.BatchNormalization()(x)
x = layers.SeparableConv2D(filters, 3, padding="same", activation='relu')(x)
x = layers.BatchNormalization()(x)
x = layers.MaxPooling2D(3, strides=2, padding="same")(x)
# DECODER
for filters in [128, 64, 32, 16]:
x = layers.Conv2DTranspose(filters, 3, padding="same", activation='relu')(x)
x = layers.BatchNormalization()(x)
x = layers.Conv2DTranspose(filters, 3, padding="same", activation='relu')(x)
x = layers.BatchNormalization()(x)
x = layers.UpSampling2D(2, interpolation="bilinear")(x)
# Add a per-pixel classification layer
outputs = layers.Conv2D(1, 3, activation="sigmoid", padding="same")(x)
# define and compile
model = Model(inputs, outputs)
model.compile(optimizer=Adam(learning_rate), loss="binary_crossentropy", metrics=BinaryAccuracy())
return model
u_net = create_unet_model(learning_rate=0.0005)
u_net.summary()
Model: "model_2"
_________________________________________________________________
Layer (type) Output Shape Param #
=================================================================
input_21 (InputLayer) [(None, 256, 256, 1)] 0
conv2d_54 (Conv2D) (None, 128, 128, 16) 160
batch_normalization_60 (Bat (None, 128, 128, 16) 64
chNormalization)
activation_4 (Activation) (None, 128, 128, 16) 0
separable_conv2d_24 (Separa (None, 128, 128, 32) 688
bleConv2D)
batch_normalization_61 (Bat (None, 128, 128, 32) 128
chNormalization)
separable_conv2d_25 (Separa (None, 128, 128, 32) 1344
bleConv2D)
batch_normalization_62 (Bat (None, 128, 128, 32) 128
chNormalization)
max_pooling2d_57 (MaxPoolin (None, 64, 64, 32) 0
g2D)
separable_conv2d_26 (Separa (None, 64, 64, 64) 2400
bleConv2D)
batch_normalization_63 (Bat (None, 64, 64, 64) 256
chNormalization)
separable_conv2d_27 (Separa (None, 64, 64, 64) 4736
bleConv2D)
batch_normalization_64 (Bat (None, 64, 64, 64) 256
chNormalization)
max_pooling2d_58 (MaxPoolin (None, 32, 32, 64) 0
g2D)
separable_conv2d_28 (Separa (None, 32, 32, 128) 8896
bleConv2D)
batch_normalization_65 (Bat (None, 32, 32, 128) 512
chNormalization)
separable_conv2d_29 (Separa (None, 32, 32, 128) 17664
bleConv2D)
batch_normalization_66 (Bat (None, 32, 32, 128) 512
chNormalization)
max_pooling2d_59 (MaxPoolin (None, 16, 16, 128) 0
g2D)
conv2d_transpose_32 (Conv2D (None, 16, 16, 128) 147584
Transpose)
batch_normalization_67 (Bat (None, 16, 16, 128) 512
chNormalization)
conv2d_transpose_33 (Conv2D (None, 16, 16, 128) 147584
Transpose)
batch_normalization_68 (Bat (None, 16, 16, 128) 512
chNormalization)
up_sampling2d_16 (UpSamplin (None, 32, 32, 128) 0
g2D)
conv2d_transpose_34 (Conv2D (None, 32, 32, 64) 73792
Transpose)
batch_normalization_69 (Bat (None, 32, 32, 64) 256
chNormalization)
conv2d_transpose_35 (Conv2D (None, 32, 32, 64) 36928
Transpose)
batch_normalization_70 (Bat (None, 32, 32, 64) 256
chNormalization)
up_sampling2d_17 (UpSamplin (None, 64, 64, 64) 0
g2D)
conv2d_transpose_36 (Conv2D (None, 64, 64, 32) 18464
Transpose)
batch_normalization_71 (Bat (None, 64, 64, 32) 128
chNormalization)
conv2d_transpose_37 (Conv2D (None, 64, 64, 32) 9248
Transpose)
batch_normalization_72 (Bat (None, 64, 64, 32) 128
chNormalization)
up_sampling2d_18 (UpSamplin (None, 128, 128, 32) 0
g2D)
conv2d_transpose_38 (Conv2D (None, 128, 128, 16) 4624
Transpose)
batch_normalization_73 (Bat (None, 128, 128, 16) 64
chNormalization)
conv2d_transpose_39 (Conv2D (None, 128, 128, 16) 2320
Transpose)
batch_normalization_74 (Bat (None, 128, 128, 16) 64
chNormalization)
up_sampling2d_19 (UpSamplin (None, 256, 256, 16) 0
g2D)
conv2d_55 (Conv2D) (None, 256, 256, 1) 145
=================================================================
Total params: 480,353
Trainable params: 478,465
Non-trainable params: 1,888
_________________________________________________________________
n_epochs = 5
model_checkpoint_2 = ModelCheckpoint(
'best_model_2.h5',
monitor='val_loss',
save_best_only=True,
restore_best_weights=True,
verbose=0
)
hist_unet = u_net.fit(train_2, validation_data=val_2, epochs=n_epochs, callbacks=[model_checkpoint_2])
Epoch 1/10 269/269 [==============================] - 19s 68ms/step - loss: 0.0225 - binary_accuracy: 0.9918 - val_loss: 0.0395 - val_binary_accuracy: 0.9860 Epoch 2/10 269/269 [==============================] - 18s 67ms/step - loss: 0.0204 - binary_accuracy: 0.9924 - val_loss: 0.0327 - val_binary_accuracy: 0.9888 Epoch 3/10 269/269 [==============================] - 18s 66ms/step - loss: 0.0179 - binary_accuracy: 0.9932 - val_loss: 0.0359 - val_binary_accuracy: 0.9874 Epoch 4/10 269/269 [==============================] - 18s 67ms/step - loss: 0.0161 - binary_accuracy: 0.9938 - val_loss: 0.0382 - val_binary_accuracy: 0.9889 Epoch 5/10 269/269 [==============================] - 18s 66ms/step - loss: 0.0151 - binary_accuracy: 0.9941 - val_loss: 0.0355 - val_binary_accuracy: 0.9885 Epoch 6/10 269/269 [==============================] - 18s 67ms/step - loss: 0.0142 - binary_accuracy: 0.9944 - val_loss: 0.0595 - val_binary_accuracy: 0.9873 Epoch 7/10 269/269 [==============================] - 18s 67ms/step - loss: 0.0132 - binary_accuracy: 0.9947 - val_loss: 0.0407 - val_binary_accuracy: 0.9877 Epoch 8/10 269/269 [==============================] - 18s 67ms/step - loss: 0.0120 - binary_accuracy: 0.9951 - val_loss: 0.0342 - val_binary_accuracy: 0.9896 Epoch 9/10 269/269 [==============================] - 18s 67ms/step - loss: 0.0114 - binary_accuracy: 0.9954 - val_loss: 0.0388 - val_binary_accuracy: 0.9875 Epoch 10/10 269/269 [==============================] - 18s 67ms/step - loss: 0.0106 - binary_accuracy: 0.9956 - val_loss: 0.0368 - val_binary_accuracy: 0.9891
4.2 - Testing the U-Net model¶
y_true = []
y_pred = []
x_true = []
for image_batch, mask_batch in test_2:
y_true_batch = mask_batch
y_pred_batch = u_net.predict(image_batch, verbose=0)
y_true.extend(y_true_batch)
y_pred.extend(y_pred_batch)
x_true.extend(image_batch)
y_true = np.array(y_true)
y_pred = np.array(y_pred)
x_true = np.array(x_true)
fig, axes = plt.subplots(nrows=3, ncols=8, figsize=(24, 9))
axes = axes.flatten()
for i in range(24):
axes[i].set_xticklabels([])
axes[i].set_yticklabels([])
for i in range(8):
axes[i].imshow(y_pred[i], cmap='gray')
axes[i+8].imshow(y_true[i], cmap='gray')
axes[i+16].imshow(x_true[i], cmap='gray')
axes[0].set_ylabel('Predicted mask', size=16)
axes[8].set_ylabel('True mask', size=16)
axes[16].set_ylabel('Image', size=16)
plt.tight_layout()
V Conclusion¶
Tumor type classification¶
We ended up with decent F1 score of 95%. Despite this fact, model still have small problem with classifying meningioma class, propably because of it's size.
Note that the dataset itself is very good quality, it's balanced in terms of specific factors like: resolution, color range, object placement.
Tumor location prediction¶
The training process takes way longer comparing to 1st model. However we can eventually increase the number of epochs to achieve better results.
Honestly I don't think that's necessary.
Significant improvement noticed after:¶
- Tweaking learning rate of the Adam optimizer
- Weighting classes
- Adding dropout layers for each convolution layer
instead of holding one dropout layer before output densing layers.
Ideas to improve:¶
- Preprocessing images with kmeans color clustering, didn't returned better results.
- Removing small images (head size),
this meant to remove images where the number of black pixels exceeded a certain %.
Overall I think this is a good starting point to create some kind of application for brain MRI scans. It'd predict tumor type and it's location.
Thanks for reading, I hope you enjoyed the analysis :]
Gracjan Pawłowski, 2024